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Abstract 

Shocks in astrophysical fluids can generate suprathermal particles by first order (or 
difl'usive) Fermi acceleration. In the test particle regime there is a simple relation 
between the spectrum of the accelerated particles and the jump conditions at the 
shock. This simple picture becomes complicated when the pressure of the accelerated 
particles becomes comparable with the pressure of the shocked fluid, so that the 
non-linear backreaction of the particles becomes non negligible and the spectrum 
is affected in a substantial way. Though only numerical simulations can provide a 
fully self-consistent approach, a more direct and easily applicable method would be 
very useful, and would allow to take into account the non-linear effects in particle 
acceleration in those cases in which they are important and still neglected. 

We present here a simple semi-analytical model that deals with these non-linear 
effects in a quantitative way. This new method, while compatible with the previous 
simplifled results, also provides a satisfactory description of the results of numerical 
simulations of shock acceleration. 
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1 Introduction 



Diffusive shock acceleration is thought to be responsible for acceleration of 
cosmic rays in several astrophysical environments. Most of the observational 
evidence for this mechanism, also known as first order Fermi acceleration, has 
been provided by studies of heliospheric shocks, but there are indirect lines 
of evidence that acceleration occurs at other shocks. A particularly impres- 
sive example was provided a few years ago by the observation of gamma ray 
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emission from the supernova remnant SN1006 [1]. These observations could 

be interpreted as inverse Compton emission of very high energy electrons, ac- 
celerated at the shock on the rim of SN1006, though other radiation processes 
may contribute [2]. 

Shock acceleration has been studied carefully and a vast literature exists on 
the topic. Some recent excellent reviews have been written [3-6]. Some of the 
more problematic aspects of the theory of particle acceleration at astrophysical 
shocks have been understood, while others are still subject of investigation. 

One of the problems that are harder to face is the problem of the injection of 
particles in the acceleration region. Only particles with a Larmor radius larger 
than the thickness of the shock are actually able to feel the discontinuity at 
the shock. The shock thickness is of the order of the Larmor radius of thermal 
protons, so that only a small fraction of the particles can be accelerated. 

The calculation of the injection efficiency is quite problematic, for several rea- 
sons: first, the distribution function of thermal particles is steeply decreasing 
with momentum, so that the number of accelerated particles changes wildly 
with changing injection momentum. Moreover the distribution of the particles 
in the shock frame is strongly anisotropic at these low momenta, which adds 
to the difficulty of obtaining a straight analytical answer. Second, the injec- 
tion of particles from the thermal distribution and their subsequent diffusive 
transport is thought to be due to the scattering against plasma waves, which 
are likely to be excited by the particles themselves, which makes the problem 
intrinsically non-linear. This non-linearity is exacerbated by the backreaction 
of the particles on the structure of the shocked fluid [7] . The only way to have 
a complete quantitative picture of the problem of shock acceleration is to use 
numerical simulations [8,6,9-11]. 

All these complicated effects, which seem to be important in several astro- 
physical situations, are nevertheless often neglected, mainly because of the 
lack of an approach that allows to take them into account without the use 
of complicated numerical simulations which are usually of restricted use. As 
a consequence, in most of the applications of the diffusive shock acceleration 
to astrophysical situations, the assumption of test particles is adopted, even 
in those cases where this approximation works poorly. Numerical simulations 
show however that even when the fraction of particles injected from the plasma 
is relatively small, the energy channelled into these few particles can be close 
to the kinetic energy of the unshocked fluid, making the test particle approach 
unsuitable. The most visible effect is on the spectrum of the accelerated par- 
ticles, which shows a peculiar flattening at the highest energies, due to the 
backreaction of accelerated particles on the fluid. The consequences on the 
spectra of secondary particles and radiation processes are clear. 
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The need to have a theoretical understanding of the non-Unear effects in parti- 
cle acceleration fueled many efforts in finding some effective though simplified 
picture of the problem. The structure of shocked fluids with a backreaction of 
accelerated particles was investigated in [12-16] in a fluid approach. The ther- 
modynamic quantities were calculated including the effects of cosmic rays, 
but the approach did not provide information on the spectral shape of the 
accelerated particles. 

In Ref. [17] a perturbative approach was adopted, in which the pressure of 
accelerated particles was treated as a small perturbation. By construction 
this method provides an answer only for weakly modified shocks. 

An alternative approach was proposed in [18-21]. This approach is based on 
the assumption that the diffusion of the particles is sufficiently energy depen- 
dent that different parts of the fluid are affected by particles with different 
average energies. The way the calculations are carried out implies a sort of 
separate solution of the transport equation for subrelativistic and relativistic 
particles, so that the two spectra must be connected at p ~ mc a posteriori. 

Recently, in [22-24], the effects of the non-linear backreaction of accelerated 
particles on the maximum achievable energy were investigated, together with 
the effects of geometry. The solution of the transport equation was written 
in [24] in an implicit form, and then expanded in terms of the unperturbed 
(linear) solution. 

Recently, some analytical solutions were also presented for the non-linear shock 
acceleration, in the particular case of Bohm diffusion coefficient [26,27]. 

The need for a practical solution of the acceleration problem in the non-linear 
regime was recognized in [28], where a simple analytical approximation of the 
non-linear spectra was presented. In this model the spectrum of the acceler- 
ated particles was assumed to consist of a broken power law, with three slopes 
characterizing the low, intermediate and high energy regimes. The basic fea- 
tures of the spectra derived from numerical simulations were reproduced with 
this method. 

In the present paper we propose an approach that puts together some of 
the elements introduced in [18], [24] and [28] and provides a semi-analytical 
solution for the spectrum of accelerated particles and for the structure of 
the shocked fluid. The method proposed is of simple use, can be adapted to 
several situations and provides results in very good agreement with numerical 
simulations, and with simplified models as that in [28]. 

The paper is structured as follows: in §2 we describe the general problem 
of linear and non-linear shock acceleration; in §3 we explain in detail our 
approach to non-linear effects in shock acceleration. In §4 we discuss the results 
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of our model and compare them with the predictions of the model in Ref. 
[28] and with the results of some numerical simulations. Our discussion and 
conclusions are presented in §5. 



2 Shock acceleration: linear and non-lineEir 



In this section we discuss the basic elements of shock acceleration and in- 
troduce our approach to the description of the non-linear effects due to the 
backreaction of the accelerated particles on the shocked fluid. 

For simplicity we limit ourselves to the case of one- dimensional shocks, but 
the introduction of different geometrical effects is relatively simple, and in fact 
many of our conclusions will be not affected by geometry. The non- linear effects 
are restricted to the mutual action of the shocked fluid and the accelerated 
particles. In other words, the present work does not include self consistently the 
production and the absorption of plasma waves by the accelerated particles. 
This simphfication is common to the approaches in [24] and [28] . 

The equation that describes the diffusive transport of particles in one dimen- 
sion is 



dt dx 



+ U— + -—P— + Q (11 
dx 3dx dp 



where f(x,p) is the distribution function, u is the fluid velocity and D is 
the diffusion coefficient. The injection of particles is assumed to occur only 
immediately upstream of the shock, so we write the source function as Q = 
Qo{p)6{x), where x — corresponds to the position of the shock front. For 
monoenergetic injection, the function Qo{p) has the following form: 



where pinj is the injection momentum and Ui is the fluid velocity immediately 
upstream {ui — u{0'^)). Ninj is the number density of particles injected at the 
shock, parametrized here as N^nj — rjNgas,!, where Ngas,i is the gas density at 
X — 0~^. The boundary condition at the shock reads: 

-P^-[D^] -[D^] +Qo{p), (3) 



3 dp \ dx 
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where U2 is the fluid velocity downstream. Here we called /o = f{0,p) the 
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distribution function at the shock position. 

A useful way of handhng eq. (1) was suggested in Ref. [24] (a similar approach 
was also adopted in Ref. [25]), and consists of integrating this equation in the 
variable x from a; = O"*" (upstream) to a; = +oo (far upstream). After some 
simple algebraic steps, in which we make use of eq. (3), we obtain the following 
equation: 



where the stationarity assumption was adopted and we assumed {Ddf/dx)2 — 
. We have introduced the quantity: 



and again we called Ui and U2 the fluid velocities at a; = 0"^ and x = 0~ 
respectively. With this formalism the compression factor at the shock is Rsub ~ 
Ui/u2- The function Up, at each momentum p has the meaning of average fluid 
velocity felt by a particle with momentum p while diffusing upstream. Since 
the diffusion is in general p-dcpcndcnt, particles with different energies will 
feel a different compression coefficient, and the correspondent local slope of 
the spectrum will be ]9-dependent. Note that, according to eq. (5), the velocity 
Up must be a monotonically increasing function of p. 

The function Up describes the mutual interaction between the accelerated par- 
ticles and the fluid. In other words, if we find the way of determining the 
function Up, as we show later, we also determine the spectrum of the acceler- 
ated particles. 

Eq. (5) is clearly an implicit definition of Up, meaning that Up depends on the 
unknown function /. However, eq. (5) allows us to extract important physical 
information. For the monoenergetic injection in eq. (2) the solution can be 
implicitly written as 









where we put (for 7^ = 5/3): 




(7) 
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Eq. (6) tells us that the spectrum of accelerated particles has a local slope 
given by 



Up — U2 Up — Imp 



The problem of determining the spectrum of accelerated particles would then 
be solved if the relation between Up and p is found This is the scope of the 
next section. 



3 The gas dynamics of modified shocks 

The velocity, density and thermodynamic properties of the fluid can be deter- 
mined by the usual conservation equations, including now the pressure of the 
accelerated particles. We write these equations between a point far upstream 
(x = +00), where the fluid velocity is uq and the density is po = fnNgas,o, and 
a generic point where the fluid upstream velocity is Up (density pp). The index 
p will denote quantities measured at the point where the fluid velocity is Up. 
We call this generic point Xp. 

The mass conservation implies: 

PqUq = PpUp. (9) 

Conservation of momentum reads: 

Poul + Pg^o = Ppul + Pg,p + PcR,p, (10) 

where Pg^ and Pg^i are the gas pressures at the point x = +00 and x = Xp 
respectively, and PcR,p is the pressure in accelerated particles at the point Xp 
(we used the symbol CR to mean cosmic rays, to be interpreted here in a 
very broad sense). In writing eqs. (9) and (10) we implicitly assumed that the 
average velocity Up as defined in eq. (5) coincides with the fluid velocity at the 
point where the particles with momentum p turn around to recross the shock. 

Our basic assumption, similar to that used in [18], is that the diffusion is p- 
dependent and that therefore particles with larger momenta move farther away 
from the shock than lower momentum particles. This assumption is expected 

^ In order to determine the spatial distribution of fluid velocity it is necessary to 
specify the exact diffusion coefficient as function of x and p. In this paper we do 
need this information, therefore the choice of D does not affect our conclusions. 
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to describe what actually happens in the case of diffusion dependent on p. As 
a consequence, at each fixed Xp only particles with momentum larger than p 
are able to affect the fiuid. Strictly speaking the validity of the assumption 
depends on how strongly the diffusion coefficient depends on the momentum p, 
but the results should not be critically affected by this assumption. Moreover, 
in case of strong shocks there are arguments that suggest that the strong 
turbulence excited by the shock should prodiicc a Bohm diffusion coefficient, 
so that the dependence D{p) on p should be at least linear. 

According to this assumption, only particles with momentum ^ p can reach 
the point x — Xp, therefore 



where v(p) is the velocity of particles whose momentum is p, and Pmax is the 
maximum momentum achievable in the specific situation under investigation. 
In realistic cases, Pmax is determined from geometry or from the duration of the 
shocked phase, or from the comparison between the time scales of acceleration 
and losses. Here we consider it as a parameter to be fixed a priori. Prom 
eq. (10) we can see that there is a maximum distance, corresponding to the 
propagation of particles with momentum Pmax such that at larger distances 
the fiuid is unaffected by the accelerated particles and Up — uq. 

We will show later that for strongly modified shocks the integral in eq. (11) is 
dominated by the region p ~ Pmax- This improves even more the validity of our 
approximation PcR.p = Pcr{> p)- This also suggests that different choices for 
the diffusion coefficient D{p) may affect the value of Pmax, but at fixed Pmax 
the spectra of the accelerated particles should not be appreciably changed. 

Assuming an adiabatic compression of the gas in the upstream region, we can 
write 



where we used the conservation of mass, eq. (9). The gas pressure far upstream 
is Pg^o = PqUq/ {'^gM'^), where 7^ is the ratio of specific heats (7^, = 5/3 for 
an ideal gas) and Mq is the fiuid Mach number far upstream. Note that eq. 
(12) cannot be applied at the shock jump, where the adiabaticity condition is 
clearly violated. 

We can rewrite eq. (10) in a convenient way, by dividing it by pqUq and using 




(11) 



p 




(12) 
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the mass conservation. We then obtain: 



1 + 



3 po^o 



Pmax 

j dp'p'Mp')fo{p') 



Pmax 



Poulpfnj 



dp'p'^v{p') exp 



Pinj 



dp" 



3upi' 



+ 



du„" 



Up" — u-2 Up" — U2 dlnp" 



where we have exphcitely written the distribution function /(p) as a power 
law with local slope Q{U), and we put U = Up/uQ. 



Differentiating the previous equation with respect to p we obtain: 

3 ( P 



dU 
dp 



I L[/-(7s+l) 



Mi 



INinjQs j P 



3 PqUq \Pinj , 



v{p) exp 



Pinj 



dp 
p 



3m. 



1 dUr. 



Up — U2 Up — U2d Inp 



or 



dU 
dlnp 



1 ^injQ. 



3 po«o 



'-vip)Pi„j exp < 



4 In 



1 - 



P 

\Pinj , 



Ml 



dp 
p 



Pinj 



SUr 



+ 



dUr 



Up — U2 Up — U2 dlnp 



.(15) 



Note that the velocity Up changes as a consequence of the pressure added by 
non-thermal particles, therefore the function U{p) must be a monotonically 
increasing function of the particle momentum. Since U{p) and dU/dlnp are 
always non-zero, we can calculate their logarithm, so that eq. (15) becomes: 



In PC/ + In 



1 L[/-(7<,+l) 



Ml 



In 



1 ^injQs I \ 

3 p^ui 



tnj 



41n 



P 

^Pinj , 



dp 

P 



^RtotU 



+ 



ot 



RtotU — 1 RtnfU — 1 



-vu 



tot^ 



,(16) 



where we put VU = dU/dlnp. The following equality is easily demonstrated: 
p 



[ 'ft 
J p 

Pinj 



dp Rtot ^ fRtotU - 1 



RtotU - 1 



R 



■sub 



(17) 



so that the equation for VU becomes: 



In VU + In 



1 ^_/-/-(79+l) 



In 



1 ^injQs I \ 

3 po^o 



+ 41n 



V 

\Pino , 



In 



RtoiV - 1 



R 



sub 



Pinj 



dp 3RtotU 

P RtotU - 1 



(18) 



Solving this differential equation provides U (p) and therefore the spectrum of 
accelerated particles, through eq. (6). 

The operative procedure for the calculation of the spectrum of accelerated 
particles is simple: we fix the boundary condition at p = pinj such that 
U{pinj) = ui/uq for some value of ui (fluid velocity at x = 0"*"). The evo- 
lution of as a function of p is determined by cq. (18). The physical solution 
must have U (pmax) = 1 because at p ^ Pmax there are no accelerated particles 
to contribute any pressure. There is a unique value of for which the fluid 
velocity at the prescribed maximum momentum Pmax is Up^^^ — Uq (or equiv- 
alently U{pmax) = !)• Finding this value of ui completely solves the problem, 
since eq. (18) provides U{p) and therefore the spectrum of accelerated par- 
ticles, calculated according to eq. (4). Conservation of energy can be easily 
checked. 

An illustration of this procedure is presented in fig. 1 where we considered the 
following special set of parameters: the Mach number far upstream is Mq = 43, 
and the gas temperature is Tq = 10^ K, corresponding to = 5 x 10® cm/s. 
The parameter 77 is taken to be 4 x 10~^, and the injection of particles occurs 
at Pinj — 10~'^'mc, where m is the mass of the accelerated particles. In this 
paper we assume that the accelerated particles are protons. The maximum 
momentum is Pmax = lO^mc. The three curves which are plotted correspond 
to the function U{p) for Rgub = 2 (dotted line), Rgub = 1-5 (dashed line) and 
Rsub — 1-817 (sohd line). It is immediately evident that only for Rgub — 1.817 
the function U(p) reaches unity at p = Pmax- Our method then provides the 
value of Rsub and consequently the values of the other parameters. 

Fig. 1 is very useful for understanding the physical meaning of the local ve- 
locity Up. Particles with large momenta feel compression factors Up/u2 which 
are larger than those felt by low momentum particles. Large compression fac- 
tors correspond to locally flatter spectra, so that the spectrum of accelerated 
particles is expected to become flatter at large momenta. We deflne the total 
compression factor Rtot — uq/u2. Note that the compression factor at the gas 
subshock is now Rsub <S 4, the value expected for a strong shock in the test 
particle regime. 
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Fig. 1. U{p) as a function of momentum, for different values of the compression 
factor at the subshock. The fluid dynamics must be unaffected at p = Pmax which 
implies that U{pmax) = 1- This condition determines the right value of Rgub- 

4 Results 



The calculations illustrated in the previous section are here tested versus pre- 
vious models and numerical simulations. In particular wc compare our results 
with the predictions of the simple approach presented in [28], that wc briefly 
summarize. In [28] the spectrum of accelerated particles has a prescribed 
shape, made of three power laws, in the low {pinj < p < mc), intermedi- 
ate (mc < p < 10~^pmax) and high {10~'^Pmax < P < Pmax) energy regimes. 
The slope in the three regions is then calculated by requiring mass, momentum 
and energy conservation. 

We first check that our model reproduces the results in the linear regime, 
where the test particle approximation can be adopted, and then we study the 
transition to the non-linear regime. For this purpose, we consider a shock with 
Mach number Mq = 5, with a gas temperature Tq = 10^ K. We choose pinj = 
lO^^mc and Pmax = lO^mc, and we study the result for different values of i], 
as plotted in Fig. 2. The sohd hue is obtained for 77 = 10~^; our model gives 
for this case Rgub — Rtot = 3.57, which is exactly the value obtained from test 
particle theory: Rsub = (8Afo/3)/ ((2/3)A'/o + 2), for 7^ = 5/3. The test particle 
approximation would provide the same compression factor for any value of 77. 
Our results for rj = 10~^ and rj = 10~^ are plotted in Fig. 2 as dotted and 
dashed hues respectively. These two cases result in {Rgub, Rtot) = (3.45, 3.83) 
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for T] = 10-3 and (Rsub,Rtot) = (2.85,4.49) for t] = IQ-^. The transition from 
unmodified shocks to strongly modified shocks is evident. 
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Fig. 2. Spectra of accelerated particles for low Mach number (Mq = 5). Increasing the 
value ofrj determines the transition from an unmodified shock (solid line, rj = 10~^) 
and a modified shock (dashed and dotted lines, with rj = 10"^ and rj = 10"^ respec- 
tively). 

In all the cases the common feature of the modified shocks is a steepening 

of the spectrum (in comparison with the hnear result) at low energy and a 
flattening at high energies. From the phenomenological point of view this can 
be of paramount importance since it may change the spectral features of the 
secondary radiation produced by the interactions of the accelerated particles. 



The more interesting case is that of strongly modifled shocks, where we expect 
the pressure in accelerated particles to become comparable with the kinetic 
energy of the upstream fluid. Intuitively this is mainly the case for high Mach 
number shocks, though we also show that there are exceptions. We can com- 
pare our results with those of the model in Ref . [28] . 

Let us consider the case Mq = 43, uq = 5 x 10® cm/s, pinj = 10-^mc, and let 
us study the resulting spectra for 77 = 10"^ and different values of Pmax- The 
results are plotted in Fig. 3. 

The results of [28] are plotted in the form of crosses {pmax = lO'^mc), stars 
{Pmax = lO^mc) and diamonds {pmax = lO^mc). The corresponding continuous 
lines are the results of our model. 
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Fig. 3. Comparison between the prediction of our model (lines) and those of Ref. 
[28] (symbols). The three sets of curves are obtained for p^ax = W'^mc (solid line 
and crosses), Pmax = lO^mc (dotted line and stars) and Pmax = W^mc (dashed line 
and diamonds). 

Some comments are in order: the slope of the spectrum of accelerated particles 
predicted in our model at p ~ Pmax is approximately equal to that obtained 
in [28] . This is not surprising since the model in [28] is based on a three power 
law approximation, and the slope of the spectrum at the highest energies is 
calculated using an asymptotic expression derived from eq. (6). At low values 
of Pmax our results are in very good agreement with the general features of 
the solutions in [28]. At increasingly larger values of Pmax the agreement is 
not extremely good but still reasonable, if one takes into account that the 
following three assumptions were made in [28]: 1) the velocity of the particles 
is assumed to be p/ m for p < mc and equal to c for p > mc; 2) the momentum 
at which there is a change in slope at low energy is forced to be at p = mc; 
3) the position of the point where there is the change in slope at intermediate 
energies is forced to be at p = 10^'^Pmax- In our model all the three assumptions 
are released and the spectra are smooth. This is the reason for the slightly 
different position of the dip in our spectra (p'^/(p)) when compared with those 
derived according to [28]. Similar conclusions hold for different values of the 
parameters. 

A more interesting comparison, to test the effectiveness of our model is that 
with the results of numerical simulations. Since a similar comparison was 
carried out in [28] to test that model, we consider here the same situation. 
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so that a full cross-check is possible. The case we consider is that of a shock 
with Mach number Mq = 128 and a fluid temperature of T = 10^ K. The 
results provided by a numerical simulation are plotted in Fig. 5 of Ref. [28] 
and reproduced in our Fig. 4 as a dashed line. From the simulation we can 
extract the values of some parameters: the injection momentum is pinj — 
7 X 10~^mc, while the maximum momentum is Pmax — lO^mc. Note that the 
simulation does not provide these parameters in a clear way, because there 
the thermal and non-thermal particles are treated in the same way, therefore 
some approximation is involved in deriving these numerical values. The best 
fit of the model in [28] to the results of the simulation implies Rgub = 2.68 and 
Rtot — 52 using 1] = 5 X 10~^. Using the same value of the injection efficiency 
r], our model predicts Rsuh = 2.365 and Rtot = 51.9. The spectra for the model 
in [28] and for our model are plotted in Fig. 4 as a solid light broken line, and 
as a thick solid line respectively. The improvement provided by our model is 
evident (note that we only used the values of rj and Pinj derived in [28] , though 
an even better agreement can be found by slightly changing these values, that, 
as we stressed above, are not clearly determined by the simulation). 




Fig. 4. Comparison between the predictions of our model (thick solid line) and the 
results of simulations (dashed line) and the approximation in Ref. [28] (solid light 
line). 

A basic issue raised in [28] is that of the accurate description of the slope 

of the spectrum of accelerated particles at p ~ Pmax- The importance of this 
point is due to the flat spectral shape at high momenta, which gives the main 
contribution to the total energy budget in accelerated particles, for strongly 
modified shocks. We show our prediction for the local slope in Fig. 5, as a 
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function of the momentum. The most noticeable feature of this figure is the 
slope ~ 3.5 at p ~ Pmax, very close to that predicted by numerical simulations. 
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Fig. 5. Slope of the spectrum of accelerated particles according with our model, with 
the parameters used to obtain the curves in Fig. 4- 

One comment is in order concerning the injection: in all the cases consid- 
ered above, the injection momentum and rj, the injection efficiency, are chosen 
independently. In a realistic case (for instance in the simulation) the two pa- 
rameters would be actually related to each other. For instance, if the particles 
are extracted from the thermal distribution, then one can write in a very gen- 
eral way, that the injection momentum is proportional to the sound speed 
downstream (cs,2), times the particle mass, pinj — ^mcs,2- The parameter ^ 
must be large enough that the particles can feel the thickness of the shock, 
determined either by the interaction pathlength of the particles, or by the 
Larmor radius of the thermal protons in the local magnetic field. Hence, ^ 
must be at least a few, and indeed it is usually assumed to be in the range 
^ ~ 4 — 10. If the expression pinj — ^mcs,2 is adopted, and the injection is 
assumed to occur from the thermal distribution, then the parameter t] is no 
longer free, and can be determined self-consistently. However it is wise to keep 
in mind that in this case r] would change wildly for small changes in ^, due 
to the exponential suppression in the distribution function at momenta larger 
than the thermal average. 

We continue now the analysis of the predictions of our model and the compar- 
ison with the results of Ref. [28]. An important issue is that of determining 
where the transition from unmodified to strongly modified shock occurs, as a 
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function of the parameters of the calculation. 

In Fig. 6 we plot the compression ratios (total and at the subshock) for different 
values of the maximum momentum of the accelerated particles. The solid lines 
represent the results of our model, compared with the values predicted by the 
simple approach in [28] (dashed lines). There is a good agreement between the 
two approaches. For rj = 10^^ the shock becomes strongly modified already for 
Pmax a few times larger than mc, confirming that even relatively low injection 
efficiencies cause the shock to be affected by the backreaction of the accelerated 
particles. 




10' 10^ 10-^ 10 10-^ 10" 10' 

LoaTp /nnc] 

Fig. 6. Compression factors according with our model (solid lines) and the simple 
model (dashed lines), as functions of pmax, for the values of the parameters reported 
in the figure. The upper curves represent Rtot O'^d the lower curves represent Rgub- 

In Fig. 7 the same compression factors are plotted versus the minimum (in- 
jection) momentum pi^j. Here some differences between our model and that 
in Ref. [28] are visible. The results can be interpreted as an evidence that our 
model predicts the onset of the modified shock phase, at slightly larger values 
of Pinj compared to the simple model. 

The results of a similar study are plotted in Fig. 8, where the dependence 
of Rfot and Rgub is investigated versus the injection efficiency rj. Again, our 
model predicts that the shock starts to be modified at values of r] slightly larger 
than for the simplified model. In particular for the values of the parameters 
reported in Fig. 8, we obtain that the shock starts to be modified at ^ 10~^, 
while the simple model gives a modified shock for 77 3 x 10~^. Finally in 



15 



25 




20 




o 


L / 7^ = 10"^ J 

/ p =10"^ mc 
/ ^max 


ctT' 

10 


/ Mq=50 
/ '^0^^ cnn/s_ 


5 





0.001 0.010 

Log[Pi^j/mc] 



Fig. 7. Compression factors according with our model (solid lines) and the simple 
model (dashed lines) as functions of the injection momentum. The upper curves 
represent Rtot O'fid the lower curves represent Rgub- 

Fig. 9 wc investigate the dependence of Rtot and Rsub on the Mach number 
of the fluid at inflnity (upstream). We speciahze our prediction to the case 
T] = 10~^, Pinj = 10~^mc and Pmax = lO^mc, but clearly similar plots can 
be produced for different regions of the parameter space. The results of our 
model are in very good agreement with the results of [28] . Our model however 
predicts a slightly lower value of the critical Mach number Mcr, above which 
the shock is no longer modifled, so that Rtot and Rsub settle down at the usual 
(linear) value of ~ 4. As found in [28], the relation between Rf^t and Mq is 
Rtot oc Mo'/", for Mo < M^,. 



5 Discussion and Conclusions 

We proposed a semi-analytical method to calculate the spectrum of particles 
accelerated at shocks, and the velocity profile of the shocked fiuid, including 
the non-hnear effects due to the backreaction of the non-thermal particles on 
the dynamics of the shocked fluid. 

The approach presented in this paper provides results that are in good agree- 
ment with the output of numerical simulations on shock acceleration and also 
with the results of a previous model aimed to a simple description of the non- 
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Fig. 8. Compression factors according with our model (solid lines) and the simple 
model (dashed lines) as functions of the injection efficiency rj. The upper curves 
represent Rtot O'fid the lower curves represent Rgub- 

linear effects [28]. Compared witli tlie latter, our model provides a better fit to 
numerical results, and is not based on a priori assumptions on the spectrum 
of accelerated particles in some energy ranges. 

We confirm the main results on shock acceleration, known from simple models 
and from simulations: i) the backreaction of the accelerated particles is impor- 
tant, even in those cases in which a small fraction of the particles injected at 
the shock is actually accelerated to suprathermal energies. For rj ^ 10""^ the 
shock becomes modified by the non-thermal pressure and both the velocity 
field of the shocked fluid and the spectrum of the accelerated particles arc af- 
fected. The minimum value of rj for which the shock is modified in our model 
is slightly larger than that predicted in [28]. ii) For linear shocks the maximum 
compression factor that can be achieved is 4 (for 7^ = 5/3) which corresponds 
to spectra oc This result is obtained for Mq 00. When the backreaction 
becomes relevant, the structure of the shocked fluid is changed into a smooth 
decrease of the fluid speed from uq at inflnity (effectively at some distance 
Dmax) to u-i at the position of the gas (ordinary) subshock. While at the sub- 
shock it is still true that the maximum compression factor is 4, the overall 
compression factor {Rtot) between downstream and Dmax can be arbitrarily 
large. Hi) The large total compression factors for strongly modifled shocks re- 
sult in a flattening of the spectra of accelerated particles at high energy. The 
slope there tends to ~ 3.5 (flatter than 4). The slope of the spectrum at low 
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Fig. 9. Compression factors according with our model (solid lines) and the sim- 
ple model (dashed lines) as Junctions of the Mach number Mq. The upper curves 
represent Rtot O'fid the lower curves represent Rgub- 

energy is determined by the compression factor at the subshock, and is usually 
steeper than the slope predicted in linear theory, iv) In general, increasing the 
Mach number corresponds to an increasingly more modified shock, up to a 
critical Mach number M^r- At M > M^r, the total and subshock compression 
factors both converge to the hnear value of 4 and the shock behaves as an 
ordinary strong shock (it is no longer an efficient accelerator). 

Although the calculations reported in the paper refer to the case of a plane 
(one dimensional) shock and the effects of waves (Alfven heating) have not 
been included, these changes are almost straightforward. In fact, for spherical 
shocks, the expression for the distribution function /o(p) remains unchanged 
[eq. (6)], although the definition of Up is formally different (it refiects the 
geometry). All the rest of the calculation remains unaffected and most of the 
results still hold. The case of spherical symmetry would be relevant for shocks 
related to supernova remnants, as discussed in [24] . 

The introduction of the heating due to damping of Alfven waves implies a 

change in the relation between Pg^p and PgX) but this change can be easily 
written in a way which is useful for our purposes. For instance, generalizing a 
discussion in [28], in the approximation of large Alfven Mach number {Ma ~ 
uo/va, where va is the Alfven speed, assumed constant in all the fiuid) we can 
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9,P 



9,0 




1 + (7. - 1) 



Ml 

Ma 




(19) 



Clearly the effects of Alfven heating can be neglected as long as <^ Ma- 
Introducing eq. (19) [instead of eq. (12)] in eq. (10), it is easy to derive an 
equation similar to eq. (18), that can be solved for U{p). The spectrum of 
accelerated particles is then obtained using the same procedure illustrated 
above. The basic effect of the Alfven heating is to reduce the total compression 
factors, and make the spectra at high energy slightly steeper than found in 
§4. 

We conclude by stressing that the non-linear effects discussed here appear to 
be relevant even for a small injection efficiency, and their phenomenological 
consequences can be critically important. For instance, as an application of the 
simple model presented in [28] , in Ref . [29] the spectra of secondary radio and 
gamma radiation in supernova remnants were calculated. The differences with 
respect to the results obtained by adopting the test particle approximation 
are impressive. 

In this prospective, the semi-analytical model presented here, being of easy 
use and having an immediate physical interpretation, provides the suitable 
tool to estimate the phenomenological consequences of non-linearity in shock 
acceleration, in the cases where it is unpractical to have access to numerical 
simulations. 
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